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o 

Q 

Q^ | We present an analytical calculation of periodic orbits in the homogeneous quartic oscillator poten- 

tial. Exploiting the properties of the periodic Lame functions that describe the orbits bifurcated 
from the fundamental linear orbit in the vicinity of the bifurcation points, we use perturbation 
Q ! theory to obtain their evolution away from the bifurcation points. As an application, we derive 

an analytical semiclassical trace formula for the density of states in the separable case, using a 
uniform approximation for the pitchfork bifurcations occurring there, which allows for full semi- 
classical quantization. For the non-integrable situations, we show that the uniform contribution of 
the bifurcating period-one orbits to the coarse-grained density of states competes with that of the 
shortest isolated orbits, but decreases with increasing chaoticity parameter a. 

m 
> 

^ ! 1 Introduction 

o. 

The homogeneous quartic oscillator potential V(x, y) = a x 4 + b y 4 + (a/2) x 2 y 2 has been the object 
[ of both classical, semiclassical and quantum- mechanical studies ||, ||, Q . Lakshminarayan et al || 

have investigated the fixed points in Poincare surfaces of section corresponding to orbits of period 
four and determined empirically some of their scaling properties. Due to the homogeneity of the 
potential in the coordinates, orbits at different energies are related to each other through a simple 
scaling of coordinates and momenta. We may therefore fix the energy E at an arbitrary value. 
The nonlinearity parameter that regulates the dynamics is the parameter a. The system possesses 
periodic straight-line orbits along both axes which undergo stability oscillations under variation of 
a. Infinite sequences of new periodic orbits bifurcate from each of these straight-line orbits and 
their repetitions, leading to almost completely chaotic dynamics || in the limit a — ► oo. 

In this paper we specialize to the symmetric case o = 6 = l/4in which the potential has C± v 
symmetry; this potential shall in the following be denoted as the Q4 potential: 

V Q i(x,y) = \(x* + y*) + ^x 2 y 2 . (1) 

The straight-line orbits along the x and y axes then obey identical equations of motion; we denote 
them as the A orbits. The dynamics of the Q4 potential (|l|) is invariant under the symmetry 
operation (cf H) a — > (3 — a)/(l + a), which corresponds to a rotation in the (x,y) plane about 
45 degrees and a simultaneous stretching of coordinates and time by a factor [2/(1 + a)] 1 / 4 . The 
limit a — * oo therefore is equivalent to the limit a. — * — 1. There are three values of a for which the 
potential is integrable: 1) a = 0, giving separability in x and y; 2) a = 1, which is the fixed point of 
the above symmetry operation, giving the isotropic quartic oscillator V(r) = \ r 4 with r 2 = x 2 + y 2 ; 
and 3) a = 3, giving separability after rotation about 45 degrees. 
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Figure 1: Stability discriminant trM of period-one orbits in the Q4 potential, plotted versus a. 
Shown are the curves for the primitive orbits A, B and C, and the orbits L and S bifurcating from 
orbits A and B, respectively. Subscripts indicate the Maslov indices a po appearing in (28). The 
vertical dashed lines at a = 0, 1 and 3 correspond to the integrable situations. 



In a recent paper M, the period-one and period-two orbits bifurcating from the A orbits have 
been classified completely in terms of periodic Lame functions. The motion of the primitive A orbit 
along the y axis is given analytically by 



x A {t) = , y A (t) = yo cn(y t, k) , y = (4S) 1/4 , 



1/2. 



with the period Ta = 4K/yo, where K = K(k) = F is the complete elliptic integral of the 

first kind with modulus k, and cn(z, k) is one of the Jacobi elliptic functions ||. The turning points 
are ±yo- Note that this solution does not depend on the value of a. The stability of the orbit A, 
however, does depend on a. The linearized equation of motion in the transverse x direction yields, 
after transformation to the scaled time variable z = yot, the Hill equation 



x"(z) + a [1 - sn 2 (z, k)] x(z) = . 



This is a special case of the Lame equation 



A"(z)+ h-n(n + l)K 2 sn 2 (z,n) A(z) = with k 



h 



a. 



-n(n+l) 



(3) 



(4) 



We therefore know here analytically the eigenvalues h = a n of the Lame equation, which correspond 
to the bifurcation points a n of the A orbit @. This agrees with the analytical result for its stability 
discriminant, given by the trace of the stability matrix M, which has been derived long ago by 



Yoshida p|: 



trM A = 4 cos ( ^Vl + 8a ) + 2 . 



(5) 

It is easily seen that the bifurcation condition trMA = +2 leads exactly to the values a n in (Q). 

In figure |l] we show the stability discriminant tr M for some period-one orbits of the Q4 potential 
in the interval — 1 < a < 4. The integrable situations are indicated by the vertical dashed lines. The 
upper horizontal dotted line is the bifurcation line tr M = +2. The orbits B are the two straight-line 
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orbits along the diagonals y = ±x, which are mapped onto the A orbits under the above-mentioned 
symmetry operation. Their motion is given by 

± x B (t)=y B (t)=y [2(l + a)]- 1 / 4 cn(y a t,K), y a = y Q [(l + a) /2} 1 '\ (6) 

their period is Tb = Ta [2/(1 + a)] 1//4 , and their stability discriminant trMe is found from trMA 
in (||) by replacing a — ► (3 — ot)/(l + a). The orbit C is a rotational orbit which has a discrete 
degeneracy of two because of time reversal symmetry; its solutions xc(t) and yc{t) could not be 
found analytically for arbitrary values of a (see section 3 for the case a = 0). L and S are the first 
librating orbits born from A and B, respectively, in isochronous (period-one) pitchfork bifurcations. 
Not shown are the period-one orbits bifurcating from B for a < —3/7. The results for trM for 
the C, L and S orbits were obtained numerically. The subscripts of all the orbits shown in figure [I] 
denote their Maslov indices a po used in the semiclassical trace formula (28) in section 3. 

The periodic solutions of the Lame equation @ are the Lame functions [||, Ec™(z, k) and 
Es™(z,k) which are even and odd functions of z, respectively, with m zeros in z € [0, 2K). For 
integer n they are polynomials of degree n in the Jacobi elliptic functions sn, cn and dn. Those with 
even m have the period T = 2K, those with odd m have T = 4K. For the special case k 2 = 1/2, m 
is fixed by m = int[(n + 1) /2] and there exists only one type of Lame polynomial, Ec or Es, for each 
value of n. It is therefore sufficient here to denote these polynomials by E n (z) with n = 0, 1, 2, . . . 
Their explicit expressions up to n = 15 have been given in Q. As there, we use from now on the 
short notation cn(z) = cn(z, k), sn(z) = sn(z, k) etc, keeping in mind that k 2 = 1/2. 

Nontrivial period doublings of the orbits A occur when trMA = —2 (cf the lower horizontal line 
in figure ||), which leads with (||) to the critical values a v = 2p (p + 1) + 3/8 with p = 0, 1, 2, . . . The 
corresponding solutions of (|j) with n = (4p + l)/2 are algebraic Lame functions (l2|, |lj| of period 
8K with m = (2p + l)/2. They are discussed in |?J in connection with the period-two orbits born 
at the corresponding (island-chain type) bifurcations. 

The purpose of this short paper is to demonstrate the use of the simple properties of the Lame 
functions for analytical classical and semiclassical studies involving the periodic orbits of the Q4 
potential. In section 2 we shall formulate a perturbation expansion for describing the evolution of 
the bifurcated orbits away from the bifurcation points a n , and in section 3 we apply some of its 
results to semiclassical calculations of the density of states. 



2 Perturbation expansion around bifurcation points 

As shown in 0], the transverse motion of the orbits bifurcated from A is, innnitesimally close to 
the bifurcation values a n , exactly described by the Lame functions. In 0, [L4j similar bifurcation 
cascades were investigated for Henon-Heiles type potentials. For these it was possible to determine 
the amplitude of the transverse motion of the bifurcated orbits from the conservation of the total 
energy, exploiting the asymptotic separability of these systems near the saddle energy into motions 
parallel and transverse to the bifurcating orbits. However, due to the scaling property of the Q4 
potential ([!]), a variation of the energy does not affect the stability of the A orbits and hence the 
energy conservation cannot be exploited in the same way. In order to find the evolution of the 
new orbits away from their bifurcations, we propose here a perturbative series expansion of the 
equations of motion around the bifurcation points a n , leading to successive analytical expressions 
for the corrections to the simple (lowest-order) Lame solutions. 
The exact equations of motion for the Hamiltonian 

H(p x ,p y ,x,y) = ^ (pl+pl) +V Q4 (x,y) (7) 

with the Q4 potential (|l|) are, in the Newtonian form, 

x + x (x 2 + ay 2 ) = 0, y + y (y 2 + a x 2 ) = . (8) 
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In the following we expand the solutions of these equations around an arbitrary bifurcation point 
a n into Taylor series 

x(t) = x A (t) + exi(t) + e 2 x 2 (t) + . . . , y(t) = y A {t) + ey 1 (t) + e 2 y 2 {t) + ... , (9) 

whereby the small dimensionless expansion parameter e is chosen, with support from numerical 
evidence, as 



a — a n \ . (10) 

As zero-order solution, we have taken XA(t) = and yA(t) = yo Ei(z) as given in (||) for the A orbit, 
since all bifurcated orbits are degenerate with the A orbit at the bifurcation points. Whereas this 
unperturbed solution is the same for all bifurcations, the corrections Xk(t), y k (t) with k = 1,2,... 
will depend explicitly on the value a n of the chosen bifurcation point. 

We now rewrite the exact equations of motion singling out the value of a n and replacing 
the difference a — a n by e 2 : 

9 9 9 Q Q Q 9 9 / \ 

oc + a n xy +exy +x =0, y + y +a n yx +eyx =0. (11) 

Hereby we have chosen a > a n which holds for all bifurcated orbits with n > 1. The orbit L3 
(n = 0) and the orbits Fq, P7 (n = 1/2, see 0) bifurcate towards smaller values of a, ie a < «o = 
and a < ai/ 2 = 3/8, respectively; for these cases the sign in front of e 2 in the equations ( |TT| ) must 
be reversed. 

Next we insert @ into (11) and extract the equations obtained separately at each order in e. 
This leads to a recursive sequence of linear second-order differential equations in the scaled time 
variable z: 

xl{z) + a n ca 1 (z)x k (z) = h k (z), y' k \z) + 3 cn 2 (z) y k (z) = g k (z) , k = 1,2,3,... (12) 

where the inhomogeneities on the rhs contain nonlinear combinations of x k i{t) and y k i{t) with 
k' < k. The homogeneous parts of these equations are identical to (||), (f|) and have the periodic 
Lame polynomials E n (z) as solutions. According to their general theory || 11 1 the second, linearly 
independent solutions are non-periodic for integer n; we shall in the following denote them by F n (z). 
We normalize them such that their Wronskians with the E n become unity: W{E n (z),F n (z)} = 1. 
The solutions E„ and F n are then related to each other by |Tl 



f z dz' f° dz' 

F »W = E » (2) i„ [l^T' M«) = F.W/ m?w . (13) 

In table |l| we give the solutions F n (z) and F n (z) for the lowest integer n. Hereby we have defined 

£{z) =2E{z) - z = f cn 2 (u)du. (14) 
Jo 

E(z) is the incomplete elliptic integral of second kind, related to that of the first kind F(ip, k) by 

E(z) = E{(p, k) = f V [l- k 2 sin 2 6} 1/2 d0 , z = F(tp, k) = [* [1 - k 2 sin 2 6}- 1/2 dO . (15) 
Jo Jo 

The function £{z) in ( p!^ ) is nonperiodic; its periodic part is given, with E = E(k) = £'(§, k), by 



Ixt (c-( : ) } ■ - -2E( ; ) -z-JL. z = 2E(z) -™z (10) 



The last equality above follows from a known relation H between the complete elliptic integrals K 
and E for k 2 = 1/2. 
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n 






F n (z) 








Eo(*) 


= Ec°(z) = 


1 




Fo(«) 


= 2 




1 


1 


Ex(z) 


= Ec|(z) = 


cn(z" 




FiW 


= 2 sn(z) dn(z) 


— cn(z) £(z) 


2 


3 


E 2 (z) 


= Es^(z) = 


dn(z] 


sn(z) 


F 2 (z) 


= z dn(z) sn(z) 


— cn(z) 


3 


6 


E 3 (*) 


= Esl(z) = 


cn(z] 


dn(z) sn(z) 


F 3 (*) 


= 2 - 3cn 4 (z) - 


- 3 cn(z) sn(z) dn(z) E(z) 



Table 1: The first four pairs of orthogonal solutions E n (z), F n (z) of the Lame equation (Q). We 
use the short notation cn(z) = cn(z, re), sn(z) = sn(z, re) etc; here re 2 = 1/2. 

For half-integer n - appearing at the period-doubling bifurcations of the A orbits leading to the 
algebraic Lame functions [12, pl| - the linearly independent solutions F n (z) are also periodic; the 



development given below must then be modified at some points. For simplicity, we limit ourselves 
in the following to the integer-n cases occurring at the bifurcations of the primitive A orbits. 
The general solutions of the equations ( |l2| ) are of the standard form 

x k (z) = c k E n (z) + d k F n (z) + H k (z) , y k (z) = a k E 2 (z) + b k F 2 (z) + G k {z) , (17) 

where the particular solutions of the inhomogeneous equations are given by 

H k (z) = I" h k (z')[E n {z)E n (z')-E n (z')E n (z))dz' 
Jo 

1 

and 

G k (z) = ( Z g k (z')[¥ 2 {z)E 2 (z')-Y 2 {z')E 2 {z)}dz' 



f z 1 f u 

E n (z) / du / h k (w)E n (w)dw (18) 



o 

= E 2 (z) [ Z du 1 ( U g k {w)E 2 {w)dw. (19) 
Jo [E 2 (u)f Jo 

The second parts of the above equations are useful for analytical computations. The coefficients 
O'k-, b k , c k and d k in ( |P7D are determined recursively by requiring x k (z) and y k (z) to be periodic 
and to have the same symmetries as the lowest-order solutions yo(z) = va(z) and x\(z) = ciE n (z), 
since these symmetries cannot be changed by varying a away from the bifurcations. The latter 
requirement leads immediately to a k = for all k. Requiring x k {z) to be periodic allows us for 
k > 3 to determine the constant c k _ 2 , appearing in different powers on the rhs of fll2|), and (for any 
k) the constant d k . Periodicity of y k {z) determines b k in terms of the cy with k' < k. We find that 
y(z) and x(z) are overall even and odd functions of e, respectively, so that c 2k = d 2k = b 2k+ \ = 
and hence x 2k {z) = y 2k +i{z) = for all k = 0, 1,2, . . . Most of the integrals can be done analytically. 

Up to this point, we have not respected the fact that the new bifurcated orbits develop their own 
periods which deviate from Ta when moving away from the bifurcation points a n . The solutions 
y k (z) and x k (z) outlined above have, indeed, all the same period Ta as the A orbit. A consequence of 
this is that these solutions do not conserve the total energy E as a function of a. In the conventional 
perturbation theory, these unwanted effects are avoided by expanding the frequencies (or periods) of 
the perturbed system in powers of e consistently along with the coordinates @. In our present case, 
this would have introduced more unknown parameters at each order, and the procedure of their 
determination would have become rather tedious. We have therefore used an alternative approach 
by an a posteriori rescaling of the dimensionless argument z = y$t in the above solutions. At a given 
order of the perturbation expansion, we set z = wt and determine the value of w by the expansion 
of the new periods 

(20) 



4K nn 

= T A 

w 



l + €T 1 + e 2 T 2 + ... 
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leading to 



w = y 



l + er 1 + e 2 T 2 + ... . (21) 



The coefficients can be determined by writing the total energy in terms of the series (||) as 

E = ±(x 2 + y 2 ) + ±(x 4 + y 4 ) + ^x 2 y 2 = E + eE 1 + e 2 E 2 + ... , (22) 

inserting the above solutions for Xk{wt) and yk(wt) with w given by (^j), and imposing that the 
energy E be independent of e (ie, of a), which means that E$ = E and E^ = for k > 0. 
Whereas the lowest-order equation Eq = E is trivially fulfilled for the solutions xa, yA in ©, the 
evaluated in terms of the above solutions Xk(yot), yk(yot) for A; > are, indeed, not all equal 
to zero. However, in terms of the rescaled solutions Xk(wt), yk(wt) we can impose E/- = for all 
k > successively by choosing the coefficients appropriately. An alternative way of deriving 
the same results (and thereby to check the algebra) is to consider the perturbation expansion (|22|) 
with new "variables" (constants of the motion) Ey. which specify the total energy E. The energy 
conservation equation ( |22|) then leads at each order k to a first-order differential equation which 
can be integrated analytically with respect to the perturbative solutions yk, with the help of ( |i~2|) 
for the Xk- From their periodicity conditions we obtain the Ej., and hence the Tf~ valid for all values 
of a and E. We found for all orbits investigated here that t\ = t 2 = t% = and T4 7^ 0, so that T 
varies with a like 

T(a) = T A [l + r 4 (a - a n f + ...], (23) 

whereby T4 is an energy-independent constant. Table ^ contains the values of T4 of the first eight 
stable orbits born at the pitchfork bifurcations of A. They will be discussed further in section 3.2. 
We see that the L3 orbit plays a special role; note that this orbit only exists at a < ao = 0, whereas 
all the other orbits given in the table exist only at a > a n . 



n 


a n 


O n 


r 4 (n) 








L 3 


37r 2 /16K 4 


3 


6 


R 5 


-11/1260 


5 


15 


L 7 


-475/354816 


7 


28 


R 9 


-4807/12627300 


9 


45 


Ln 


-160425/1092591808 


11 


66 


R43 


-75981/1115482060 


13 


91 


L15 


-4626964303/129311102954880 


15 


120 


R17 


-56892225 /2767492640836 



Table 2: Lowest non- vanishing perturbation term T4 in the period (23) of the first eight stable 
orbits O n born at the bifurcations of the primitive orbit A. 

As an illustration of our method, we present the results for the orbit R5 born at 03 = 6 and 
existing only for a > 6. We obtain the following analytical solutions up to order e 4 (with z = wt) 



X3(t) = — c\ sn(z) dn(z) cn(z) 
~2 



x\(t) = c\ sn(z) dn(z) cn(z) , 

y {t) = wcn(z) = yA{z), 2/2 (*) 

y 4 (t) = £icn(z){— 1— [-797- 206cn 4 (z) + lllcn 8 (z 
w I 10 080 L 



cn (z) 

132 180 v ; 



4 ^cn(z) [l+cn 4 (z)], 



+ 



264 



l + cn 4 (z) 



(24) 



where d = Wy/lT/45, and r 4 = -11/1260 so that w = y Q /[l - (ll/1260)e 4 ]. 
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Figure 2: Coordinates x{t) and y(t) of the orbit R5, evaluated at a = 7, ie e = 1 (left side) and 
a = 8, ie e = \/2 (right side). Solid lines: numerical results. Heavy dotted lines: perturbation series 
(||) including the terms in ( |24| ) up to 3rd and 4th order for x(t) and y(t), respectively. Other lines 
give partial contributions as labeled by the inserts. (Units such that h = 1; E = 1/6.) 

In figure [2] we compare the results obtained by summing the above analytical results according 
to (H) up to order e 4 , shown by the heavy dotted lines, with the results obtained by numerical 
solution of the exact equations of motion (H), shown by the solid lines. The single contributions 
obtained at each order in e are shown by the other lines, as labeled by the inserts. The left panels 
are calculated for a = 7 where e = 1, and the right panels for a = 8 where e = y/2. The convergence 
is surprisingly good even when the expansion parameter e is larger than unity. This is due to the 
rapidly decreasing amplitudes of the Xk(z) and yk(z) with increasing k. Note that for y(t), where 
the agreement is perfect for e = 1 and still very good for e = y/2, we have included three terms, 
whereas x(t) only contains two terms. Obtaining x^(z) would have required to calculate both ye(z) 
and x-j(z) and to make them periodic, which - though analytically possible - would have been 
rather cumbersome. But the fact that the remaining errors in x(t) are of the same order as e^y^z) 
suggests that adding the third term e 5 xs(z) would lead to an equally good convergence for x(t). 
Similar results were also obtained for other bifurcated orbits. 

For the period-two orbits born at period-doubling bifurcations, the algebraic Lame functions 
have to be used. The repeated integrations arising in the perturbation expansion then become more 
difficult and we could not do all of them analytically. Resorting to numerical integrations, however, 
whereby the coefficients bk, Ck and dk in (|T7|) were determined numerically by iteration, we could 
reach a similar convergence of the perturbation series. 
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3 Semiclassical trace formulae for the density of states 



In this section we shall apply our perturbative results to investigate the role of the pitchfork bi- 
furcations of the A orbit in semiclassical calculations of the density of states. We coarse- grain 
both the exact and the semiclassical density of states by a convolution with a normalized Gaussian 
exp{— (.E/7) 2 }/7\/7r. The quantum-mechanical coarse-grained density of states is then given by 



g qm (E) = J-]Texp{-[(£-£„)/ 7 ] 2 } 

7V 71 " n 



(25) 



in terms of the exact quantum spectrum {E n } which we have obtained by diagonalization of (|l]) in 
a harmonic oscillator basis. For larger values of 7, the prominent gross-shell structure in the density 
of states is emphasized while finer details of its oscillations are suppressed. In the limit 7 — > the 
sum of delta functions is recovered. The oscillating part of (|25|) is defined as 



Sg qm (E) = g qm (E) - g TF (E) , (26) 

where gTF(E) is the average part obtained in the Thomas-Fermi approximation which for the Q4 
potential can be calculated analytically: 

2 



9tf{E) 




for - 1 < a < +1 , 



for a > +1 



(27) 



Higher-order fi corrections to the average density of states |T7[] are negligible in the present system. 
The semiclassical trace formula for the density of states has the general form 

"1 

Ti 



5g sc (E) = -p^2\o(E) / 7 (T po ) cos 



po 



S P o(E) 



7T 



a 



po , 



(28) 



The sum goes over all periodic orbits (po) of the classical system. S po (E) are the action integrals 
§ p • dq along the periodic orbits and a po the so-called Maslov indices. The amplitudes A po (E) 
depend on the number of constants of motion of the system. The factor / 7 (T po ) is given by |l7]] 



po! 



exp 



(nT po /2hf 



(29) 



where T po = dS po {E) / dE are the periods of the orbits. This coarse-graining factor favours the 
contributions of the shortest orbits to the gross-shell structure obtained with larger values of 7. 
When the system possesses no other constant of the motion besides the energy, all orbits are 



isolated and Gutzwiller's original form [16| of the amplitudes A po (E) applies. In the presence of 
continuous symmetries, and for integrable systems in general, other forms must be used (see |l7j 
for a survey of trace formulae and the calculation of A po and a po ) . 

In the presence of bifurcations, the Gutzwiller amplitudes A po (E) for isolated orbits cannot be 
used. However, the uniform approximation for generic pitchfork bifurcations 19 1 can be applied 
to describe the isochronous bifurcations of the A orbit discussed in this paper. The factor 2 appearing 
in the period doubling for the generic case here plays the role of the extra degeneracy factor 2 of 
the bifurcated orbit pairs, which is due either to their reflection symmetry at the symmetry axis 
containing the A orbit (for the librating orbits L n ) or to the time reversal symmetry (for the rotating 
orbits R n ). Adapting the results of |n| to the present system, the combined contribution to the 
semiclassical density of states from all the orbits participating in the bifurcation at a n (including 
their degeneracy factors) is given by 



5g b sl f (E,a r 



r A r(i/4) 



COS 



kS A 



- (-1)"^ - sign(a n )| 



(30) 
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Here k is the repetition number of the orbits and v n the Maslov index, corresponding to a po in 



(28), of the unstable primitive orbit involved. The parameter a n stems from the normal form of 



the action function S(q',p) used in the phase-space representation of the trace integral at the n-th 



bifurcation of the A orbit. From the expansions given in [19] for the properties of the periodic orbits 



near the bifurcation, and using our results in the previous section, we can determine a n to be 



1 / 4vr \ 2 1 



° n -4S A U + J r 4 (n)' (31) 



In the following, we want to examine the importance of the contribution ( |30|) in relation to that 
of other non-bifurcating orbits of the system. Quite generally, in an integrable two-dimensional 
system the leading families of degenerate orbits have an amplitude proportional to H~ 3 ^ 2 (except, 
eg, for harmonic oscillators with rational frequency ratios, where the leading amplitudes go like 
h~ 2 , see p0[]). They are therefore expected to dominate over the bifurcating orbits which here have 
an amplitude proportional to ?i _5//4 as seen in (|30|). In section 3.1 below we will test the relative 
importance of the latter for the particular integrable situation given in our system for a = 0. On 
the other hand, the standard Gutzwiller amplitudes for isolated orbits in a non-integrable system 
go like h~ , so that the bifurcating orbits can be expected to play a larger relative role. As we will 
see in section 3.2 their relative weight is, however, subject to a rather subtle balance between the 
stability of the shortest isolated orbits and the a dependence of the factor |a n | _1//4 in the amplitude 
of ©. 

3.1 The separable case a = 

For a = the potential (jl]) is separable and thus defines an integrable system in which the majority 
of the periodic orbits live on a 2D torus in phase space. The EBK quantization can therefore be 
applied to obtain a semiclassical quantum spectrum from which a trace formula is readily derived 



[ 21 1 . However, the system possesses also the isolated orbit A which undergoes a bifurcation at a = 0. 
As we see from figure |l] the orbit A, which is stable (A3) for < a < 1, becomes unstable (A2) for 
a < 0; the orbit L3 bifurcating from it exists only for a < and is stable for —0.5 ^ a < 0. This 
gives a nontrivial contribution to the density of states whose role will be investigated numerically 
below. 

Straightforward EBK quantization of the 2D torus gives the approximate spectrum 

= (1/4) (6^/4K) 4 / 3 \(n x + 1/2) 4 / 3 + (n y + 1/2) 4 / 3 ] . (n x ,n y = 0,1,2,...) (32) 



From this we obtain in the standard way [21, 22 1 the following Berry- Tabor type trace formula to 
leading order in U: 



T Snm{E) - — 



(33) 



n,m=l 

The classical actions and periods 

S nm {E) = (4K/3) (4£) 3 / 4 (n 4 + m 4 ) 1 / 4 , T nm {E) = 4K (AE)- l ' A (n 4 + m 4 ) 1 / 4 (34) 
of the 2D rational tori (n, m) correspond to the periodic orbits 

x(t) = w x cn(w x t + cj)) , y(t) = w y cn(w y t) , (35) 

with 

Wx = (4£) 1/4 n(n 4 + m 4 )- 1/4 , w y = (4£) 1/4 m (n 4 + m 4 )" 1/4 . (36) 

These orbits form degenerate families described by the parameter <fi £ [0, 4K). The (n, n) resonances 
are the families containing the n-th repetitions of the orbits B (for <p = or 2K) and C (for <j) = K 



or 3K). Summation over all resonances (n,m) with n,m > yields the trace formula (|33] 
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However, the system also contains the isolated resonances (n, 0) and (0, m) with n,m= 1,2,3,... 
which correspond to the (repeated) one-dimensional A orbits in the x and y direction. With the 
uniform approximation ( |30| ) discussed above, we can include them in the trace formula, together 
with the orbits L3 born at the bifurcation A3 —* A2 + L3. This gives the following common 
contribution of all A and L orbits and their k-th repetitions 

(4K) 3 / 4 



7T 



00 

-1/16 E( 
k=l 



1 



cos 



S A (E) 



3tt 



(37) 



where S\(E) = (4K/3) (4I?) 3 / 4 is the action of the primitive A orbit. Hereby we have used the 
value r 4 = 3vr 2 /16K 4 for the L 3 orbit as given in table |, and the relation [| T(l/4) = (i^K) 1 / 2 . 

The total oscillating part of the semiclassical density of states at a = then is given by the 
sum of the contributions ( |33"|) and (|37|). We observe that the latter is of order Tl 1 ^ relative to the 
leading-order contributions of the 2D torus families. (For isolated orbits, this relative factor would 
be h 1 / 2 as mentioned above.) We therefore expect that the bifurcating orbits have a non-negligible 
influence on the density of states, at least for low energies where the negative power of E in the 
amplitude of (|37|) does not suppress their contribution too much. 

In figure |3] we show the coarse-grained level densities obtained with 7 = 1 (in units such that 
1); the periodic orbit sums in p3|) and (37) can here be limited to n mo 



2. In 



the upper panel, we show separately the semiclassical results (33) of the 2D tori (solid line) and 
of the bifurcating orbits (dotted line). Both are seen to have a monotonously increasing amplitude of 
the oscillations, whereas the quantum result (dashed line) exhibits a pronounced beating structure. 
When adding both semiclassical contributions, as shown in the lower panel, the quantum result is 
nicely reproduced. The quantum beat is clearly the result of the interference between the shortest 
orbits of the torus family on one hand and the bifurcating isolated orbits A and L3 on the other 
hand. 
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Figure 3: Oscillating part of the coarse-grained density of states at a = (units such that 
h = 1; Gaussian smoothing width 7 = 1). Dashed lines: quantum-mechanical results (26). Solid 
and dotted lines: semiclassical results with n max = m max = k max = 2. Upper panel: separate 
contributions ( |33"| ) of 2D tori and ( |37[ ) of bifurcating orbits A + L3. Lower panel: the solid line is 
here the sum oiDoth semiclassical contributions. 
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Figure 4: The same as figure ^, but for the total densities of states including the average part (27), 
with 7 = 0.06. Periodic orbit Slims 3X6 Cllt clt Tl max — 771 
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30. 



The effect of the bifurcating orbits is less dramatic when we choose a finer resolution of the 
energy spectrum, obtained with a smaller value of 7. In figure |] we show the total density of 
state including the average part ([27]), obtained for 7 = 0.06. Hereby periodic orbits up to n max = 
f^max = k max ~ 30 contribute to the semiclassical results. We have normalized g {E) here by the 
factor 7y / 7T so that the quantum result exhibits the correct degeneracies 1 or 2. Some apparently 
wrong higher degeneracies are badly resolved accidental (near-)degeneracies. At first glance there is 
not much of a difference between the upper and lower panels of figure However, a closer inspec- 
tion reveals that the semiclassical degeneracies are wrong, with errors of up to ~ 20%, in the upper 
part where only the 2D tori are included. They come much closer to the exact ones when including 
the bifurcating orbits, as seen in the lower panel. There also the bottom regions between the peaks 
are improved, the remaining small oscillations being numerical noise. Note that the inclusion of 
the bifurcating orbits does not affect the semiclassical peak positions which are exactly those of the 
spectrum -E™ K given in (^). The shifts in the peaks seen at the lowest energies are due to the 
typical errors inherent in the EBK approximation, which rapidly decrease with increasing energy. 



3.2 The non-integrable case a = 6 

As an example of a non-integrable situation we choose a = 6 where the rotating orbit R5 is born 
at the bifurcation A5 — > A$ + R5. In figure || we show a Poincare surface of section (y,p y ) taken 
at x = 0. The large regular island in the middle contains the orbits A and R5 at the origin. At 
its border we see an island chain of eight pairs of stable and unstable orbits born from a period- 
quadrupling bifurcation of the A orbit at a = 5.4306; the scaling properties of these fixed points 
were investigated in ||. The four drop-like small regular islands around (y,p y ) = (±0.431, ±0.362) 
contain the stable librating orbits S3 that are born from the B orbits in the isochronous pitchfork 
bifurcation B3 — > B2 + S3 at a = 3, as seen in figure [[]. [This is the same scenario as A3 —* A2 + 
L3 at q = 0, obtained after the transformation a — > (3 — a)/(l + a) discussed in the introduction.] 
The orbits S3 are still stable at a = 6. Since their common periods and actions are smaller than 
those of the orbit A, we expect them to influence the gross-shell structure in the density of states. 
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Figure 5: Poincare surface of section (taken at x = 0) for a = 6. See text for a discussion of the 
visible regular islands. 

Clearly, there is no hope of obtaining full semiclassical quantization by summing over all periodic 
orbits in this mixed phase-space system, where many orbits of longer periods will also be close to 
bifurcations. However, we may try to obtain the coarse-grained density of states using only the 
shortest periodic orbits, as in figure ||| above, and thus using a cut-off k max of the sum over the 
repetition number k. The contribution of the four S3 orbits to the semiclassical trace formula has 
the standard form for a stable isolated orbit fig, B3]: 



= - h T S3 (E) £ /#T S3 ) S in 



(38) 



whereby the period Ts 3 (E), the action integral Ss 3 (E) and the stability angle xs 3 have been de- 
termined numerically. The contribution from the bifurcating A and R orbits is given by (|30| ) with 
v n = 6 and using J3l| ) with T4 = —11/1260 according to table |2| for n = 3. 

Figure ^ shows the coarse-grained density of states obtained with a Gaussian smoothing width 
of 7 = 1.0, including only the k max = 2 lowest harmonics of the shortest periodic orbits. The upper 
panel gives the separate contributions of the isolated S3 orbits (solid line) and of the bifurcating A 
and R orbits (dotted line); the lower panel gives their sum. The dashed lines give the quantum- 
mechanical result in both panels. As in the case shown in figure pi each of the single contributions 
give a monotonically increasing amplitude of the oscillations in 5g sc (E), and only their superposition 
yields a beating result that reproduces the quantum-mechanical result fairly well. The remaining 
discrepancies must be due to the missing contributions from other orbits with longer periods, of 
which there are already many at a = 6. 

In spite of the fact that their amplitude in ( |30|) is larger by a relative factor 7l _1//4 than that of 
the isolated orbits, the overall contribution of the bifurcating orbits to the coarse-grained density 
of states is seen in the upper part of figure ^ to be lower. This has two reasons. First, the shortest 
isolated orbits S3 have a period that is substantially smaller than that of the orbits A and R5 , which 
leads to a larger value of the coarse-graining factor / 7 in the trace formula. Second, the orbits S3 
are stable and have therefore a stability denominator of order unity, whereas the amplitude in (|30|) 
has a larger overall denominator. Still, the uniform A + R5 contribution is strong enough to cause 
the beating interference pattern seen in the total 5g sc (E). 
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Figure 6: Oscillating part of density of states as in figure |, but at a = 6. Upper panel: separate 
contributions of the stable isolated orbit S3 (solid line) and of the bifurcating orbits A + R5 (dotted 
line). Lower panel: sum of orbits S3 + A + R5 (solid line). Truncation of the trace formulae (|30|) 
and ( p8| ) at k max = 2, Gaussian averaging width 7 = 1.0. The dashed line in both panels gives the 
quantum-mechanical result. 

It is now interesting to speculate about the relative contribution of the bifurcating orbits when 
the parameter a is further increased on the route towards chaos. On one hand, the average stability 
of the shortest isolated orbits can be expected to decrease. On the other hand, figure ^ reveals that 
the uniform amplitude in (p0|) of the period-one (k = 1) orbit A with its offsprings at the bifurcation 
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Figure 7: Amplitude of the uniform approximation (^) for the bifurcating A orbit (k = 1) to the 
semiclassical density of states (E = 1/6). The crosses give the values obtained using the values of 
T4 in table for n > at the bifurcation points a n ; the solid line is an asymptotic fit for large a. 
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points also decreases with increasing a n . (The crosses give the analytical values obtained using 
the T4 values from table [2|, and the solid line is a numerical fit showing that the amplitude decays 
asymptotically like a _1//4 .) Taking furthermore account of the fact that more and more new pairs of 
isolated orbits are created from the bifurcations along the route to chaos (and new orbits bifurcate 
from these again), even the number of shortest (period-one) isolated orbits increases very fast. We 
therefore tend to conclude that the relative weight of the bifurcating orbits in the coarse-grained 
density of states will decrease with a. It would, however, require a major numerical effort to study 
this balance in more detail and, in particular, to investigate the contributions from higher repetitions 
to the finer details of the density of states. 



4 Summary 



We have formulated a perturbative scheme to calculate analytically the shapes of periodic orbits 
bifurcating from the straight-line orbits in the quartic oscillator potential, exploiting the nice prop- 
erties of the periodic Lame functions in terms of Jacobi elliptic functions. In order to simplify 
the recursive determination of the unknown constants in the perturbation expansion, we have used 
an unconventional a posteriori rescaling of the perturbed solutions that takes care of the varying 
periods of the bifurcated orbits. 

For stable period-one orbits we are able to give analytical expressions of the resulting perturba- 
tive series. Even when the perturbation parameter e is greater than unity, satisfactory convergence 
to the numerically obtained solutions can be reached by going up to order e 4 where the necessary 
algebraic work is still not too demanding. 

We have studied semiclassical trace formulae for approximating the quantum-mechanical density 
of states of the quartic oscillator. In the separable case a = we could give a completely analytical 



trace formula, consisting of the contributions (27), (p3|) and (|37|). It was shown numerically to give a 
very good approximation to the exact quantum-mechanical density of states. Its low-frequency part, 
extracted with an energy coarse-graining parameter 7 = 1, reveals a strong quantum beat. This 
beat can only be reproduced semiclassically when the shortest periodic orbit families are allowed to 
interfere with the isolated orbits A and L3 taking part in a pitchfork bifurcation. The contribution 
of the latter was calculated using a slightly modified uniform approximation for generic pitchfork 



bifurcations given in (37). The high-resolution spectrum is dominated by the Berry- Tabor part 
which yields exactly the peaks corresponding to the EBK spectrum. The quantum degeneracies 
are, however, substantially affected by the contributions from the bifurcating orbits. 

In the non-integrable situation a = 6, which exhibits a strongly mixed phase space, we could 
also reproduce the coarse-grained quantum density of states semiclassically, including besides the 
bifurcating period-one orbits A and R5 also the shortest stable isolated orbits S3 which had to be 
obtained numerically. The two types of orbits interfere, again leading to a beat structure. Finally, 
we have shown the uniform contribution of the bifurcating A orbits and its offsprings to decrease 
like a -1 / 4 for further increasing values of the chaoticity parameter a. 

With these results we have demonstrated the possibility of accessing analytically the periodic 
orbits in a system with mixed classical dynamics. We hope to stimulate further research in this 
direction, in particular towards semiclassical studies of the role of periodic orbits and their bifur- 



cations in connection with spectral statistics |y, 24|, which would have exceeded the scope of our 



present short report. We also aim at further improving our understanding of the uniformization 
of semiclassical trace formulae, including bifurcations of higher codimensions and - as a long-term 
(and perhaps too ambitious?) goal - including complete bifurcation cascades. 
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Note added in proof 



We have confirmed the result (^) m an alternative way, without making use of any normal form. 
We start from the Poisson summation-integration of the density of states Y^n x ,n y &{E — E^^) using 
the EBK spectrum ([^), like in the standard derivation [21] of (^). The first integration can be 
done exactly exploiting the delta function. For the edge contributions corresponding to the A and 
L3 orbits, the usual stationary phase approximation cannot be used and the second integration must 
be done more carefully. Its asymptotic evaluation, valid here for kSx/Ti > 2, leads precisely to the 
analytical result (37). 
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